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A simple approximate solution to the linear response equations of time- dependent density func- 
tional theory (TDDFT) is given. This extends the single-pole approximation (SPA) to two strongly- 
coupled poles. The analysis provides both an illustration of how TDDFT works when strong 
exchange-correlation effects are present and insight into such corrections. For example, interac- 
tion can cause a transition to vanish entirely from the optical spectrum. 
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INTRODUCTION 

Ground-state density functional theory (DFT) has 
been very successful for atoms, molecules, and solids 
P, Q- Similar success is now being enjoyed by time- 
dependent DFT (TDDFT) H 0, because of its com- 
bination of accuracy combined with low computational 
cost |5|. While TDDFT has a huge variety of applica- 
tions 6], it is low- lying photo-excitations of molecules 
that has seen its greatest use ||. 

In the present work, we restrict our discussion to lin- 
ear response of a non-degenerate ground-state. Just as in 
ground-state DFT, all many-body effects, i.e., exchange 
and correlation (XC), are contained in a well-defined 
functional, the XC kernel Q. In any practical calcula- 
tion, this functional must be approximated. In most cal- 
culations, an adiabatic approximation is made, and the 
static limit of the kernel is applied. Typical approxima- 
tions are then adiabatic local density approximation Q 
(ALDA) or generalized gradient approximation, or exact 
exchange U0> [ZZ3| • The reliability and accuracy of these 
approximations to TDDFT is much less well-understood 
than it is in ground-state DFT. 

One can do many calculations on many systems, in 
order to gain insight into the accuracy and reliability 
of theory, but it can be much more effective to develop 
simple approximations to the solution of the TDDFT re- 
sponse problem [111 . A classic example is the single-pole 
approximation j!2j . within which TDDFT yields a sim- 
ple correction to the KS transition frequencies which is 
just the expectation value of the Hartree-XC kernel on 
the transition orbitals. While usually accurate 0], the 
most important feature of this approximation is the in- 
sight it yields into the workings of TDDFT. It yields a 
first approximation to TDDFT effects with almost no ex- 
tra effort beyond a ground-state calculation, and gives a 
simple picture for such effects 0] . It has also been shown 
[]t1 | that, if a transition is only weakly-coupled to others 
in the system, one can use this to estimate the XC kernel 



itself. Unfortunately, this is rarely the case in practice. 

In the present work, we generalize the SPA to a double- 
pole approximation (DPA), in which we explicitly solve 
the TDDFT response equations for exactly two transi- 
tions. This produces a variety of results beyond that 
of SPA. Most importantly, one can study TDDFT XC 
corrections to KS levels when there is strong coupling 
between levels. But one can also see when SPA fails, and 
recover Gorling-Levy perturbation theory |14| results for 
the coupling-constant expansion of excited states [Tlj . 
DPA has recently been successfully applied to core-hole 
interaction in the X-ray absorption spectroscopy of 3d 
transition metals [lfij |. 



DOUBLE-POLE APPROXIMATION 

In the matrix-formulation of the TDDFT response 
equation within the adiabatic approximation, the exact 
eigenvalues and oscillator strengths can be obtained from 
the solution of the following eigenvalue problem 0] 

^VlV(f!) iy =tt 2 v q , (1) 

q' 

where the matrix W is given by 

Wgq, (ft) = S qq , + 4 M m , (SI) (2) 

and 

M qg >(Sl) = / d 3 r J dV$*(r)Jf(rr'fi)$ g ,(r'). (3) 

Here iOi is the Kohn-Sham transition frequency and for 
single particle transitions q (q = k — > j) the shorthand 
$ 9 (r) := <pfc(r)<^*(r) has been introduced. The kernel 
K(r, r',oj) consists of the bare Coulomb interaction and 
the approximate XC kernel / xc (r, r', w): 



K(t,t',u) = - + / xc (r, r', u). 

\r — r'| 



(4) 
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Atomic units (e 2 = h = m = 1) are used throughout. 

We now solve these equations exactly for a 2x 2-system, 
i.e., ignoring coupling to all other transitions. To simplify 
the discussion we assume a frequency independent kernel 
and real orbitals, i.e. M qq i = M q / q . Thus the relation 
between matrix elements of Casida's equation and the 
kernel is: 

W u =w? +4wjAf«, W 12 =4V^I^Mi 2 . (5) 
Next define the average 



W = 1 (W u + W 22 ) 



(6) 



and difference 

AW = W22 - Wu (J) 
of the diagonal elements. We define a mixing angle by: 

2W 12 



tan 8 = 



AW 



(8) 



choosing the branch between and ir. The eigenvalues 
can then be written succinctly as 

_ 1 AW 



± 2 cos6» 

while the normalized eigenvectors are 
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(9) 



(10) 



The physical oscillator strength can be obtained from the 
following expression 0] 



where 




x 2 



(11) 



(12) 



and the x^ s denote dipole matrix elements of KS or- 
bitals. Given that there are only two transitions, we give 
a geometric meaning to the oscillator strengths. Writing 

A KS = sin 2 a KS , / 2 KS = cos 2 « KS (13) 



and 



we find 



/_ = sin 2 a, /+ = cos 2 a, 



a = a KS - 6/2, 



(14) 



(15) 



i.e., the oscillator strengths are represented by a unit 
vector in 2d space, and the coupling merely rotates this 
vector. Note that the Thomas-Reiche-Kuhn (TRK) sum 
rule (sum of the oscillator strengths is 1) is obviously 
preserved. 



SINGLE-POLE APPROXIMATION 

As mentioned above, the single-pole approximation is 
a useful approximation to TDDFT results. We recover 
SPA results by inserting 9 — in our formulas. Thus 



SPA 



— AW 
W± — 



(16) 



and the oscillator strengths reduce to their KS values. 

We can now study the leading corrections to SPA pro- 
duced by DPA when the coupling between poles is weak. 
Writing 77 = W12/ 'AW , and assuming rj << 1, for the 
eigenvalues, we find 

n± = nf A ±^ xV + o( v % (17) 



while for the oscillator strengths, we have 

/+ = /f s + 2, ?V / /P/P + 0( 7? 2 ), 

/_ = fr-^\/f2 Ks fi KS +o(v 2 ). 



(18) 



Note that the corrections to the peak positions are 
second-order in W12, while the corrections to \[J± are 
first-order. Thus SPA is expected to be much better for 
peak positions than for peak heights. 

Lastly, we point out that this expansion was deduced 
for the general case in Ref. 01 j an d used (among other 
things) to identify coefficients in the Gorling-Levy ex- 
pansion of excited-state energies. Our results here agree 
with those, but in the special case of transitions to which 
DPA applies, yield results that include a resummation of 
all orders in the adiabatic coupling constant of density 
functional theory. 



HIGH-FREQUENCY LIMIT 

So far we have given exact results for the double-pole 
approximation. However, in many cases where DPA ap- 
plies, there is a further simplification. Usually the two 
transitions are closer to each other than any others that 
couple to the pair. If in addition their frequency differ- 
ence is small relative to their mean frequency, for both 
the interacting and KS systems, i.e., 



Q, uj » Afi, Alu, 



(19) 



we find much simpler results, which are very useful for 
interpretation. 

The SPA discussed above reduces to 



nf A = io t + 2M u , (i = 1,2). 



(20) 



In fact, the original SPA was applied for just a forward 
transition, yielding exactly this result 12]. However, this 
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FIG. 1: Interacting and Kohn-Sham spectra as function of FIG. 2: The scaled coupling angle 8/tt as function of the 



frequency (wi = 9 eV, M12 = 0.2 eV). 



position of the lower transition. 



approximation violates the TRK sum rule when applied 
outside the high-frequency limit, so the symmetric result 
(sometimes called the small-matrix approximation [l?! 
llSj ) is preferable. We use the term SPA to mean the 
symmetric result throughout this paper. 
The mixing angle is given by 



tan# = 



4Mi 2 
Afl SPA ' 



(21) 



i.e., it is the ratio of the off-diagonal matrix elements 
of the kernel on the scale of the separation in SPA that 
matters. We find 



r> SPA ± 



2cos6> ' 



(22) 



SPA yields the correct average position of the two lines, 
but their splitting is greater than SPA predicts (level 
repulsion). 



In the rest of this section, we explore what happens in 
the DPA model of TDDFT. In order to emphasize that it 
is not the absolute magnitude of the off-diagonal matrix 
element that is significant, but rather its strength relative 
to the separation between the peaks, we now consider 
all the same parameters, but imagine increasing 10 \. In 
Fig. [21 we plot the mixing angle as a function of oj\. 
At uj c = 2 (-3 + V69) eV « 10.61 eV, the diagonal 
matrix elements Wu match, so that AW — and 9 = 
tt/2. At that point, the peaks are a 50:50 mixture of 
the two KS levels. In that region, the levels are strongly 
coupled, and the spectrum distorts mightily from its KS 
shape. The width of the transition region can be defined 
as the change in frequency needed to bring 9 from 7r/4 
to 37r/4, and, from Eq. H21(l in the high-frequency limit, 
is seen to be 4Mi2 = 0.8 eV here, i.e., proportional to 
the off-diagonal element, but quite a bit larger. More 
significantly, there are tails in the transition that decay 



ILLUSTRATIONS 



To illustrate our results, consider a weak lower- 
frequency transition (u>i — 9 eV, f^ s = 1/10) and a 
strong higher-frequency transition {u)2 — 12 eV, f 2 S = 
9/10). We imagine these have significant diagonal kernel 
matrix elements M%i = 3 eV, M22 = 2 eV, but are not 
strongly coupled to one another, M 12 = 0.2eV. We have 
plotted the interacting and KS spectra in Fig. ^ The 
peaks are Lorentzians of width 0.2, mimicking a mea- 
surement of finite resolution. Because the coupling is 
weak, the single-pole approximation is excellent, and ac- 
curately predicts the large shifts in positions. However, 
SPA wrongly predicts no variation in oscillator strength. 
In fact, one can see from the figure that the first peak 
has actually lost intensity relative to its KS value. 
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FIG. 3: Interacting and Kohn-Sham excitation energies as 
function of u>\. 



extremely slowly with pole separation. On the contrary, 
SPA yields a function that steps from to 1 at 11 eV. 

To see this, in Fig. |3J we plot the interacting levels 
n± as a function of u>i, and observe the avoided cross- 
ing. Note that straight line plots, extrapolated from the 
limits where lo\ is either far above or far below w c , yield 
extremely accurate results almost everywhere. This is 
the SPA result. In fact, from Eqs. 10 and ((TBI) , we see 
that the crossover point is exactly given by SPA. More- 
over, in the high-frequency limit, Eqns. (|20I22|) yield 



1.0 



I An 



SPA 1 2 



16|Mi2 



(23) 



So if the off-diagonal matrix elements are small relative 
to the SPA separation, the true separation is not much 
greater; the closest the two levels come is a separation of 
4|Afi2|, i.e., they never cross. 

But in Fig. 01 we plot the associated oscillator 
strengths. The effect of coupling is extremely dramatic. 
Note first that, for oj\ below the strong coupling region, 
the bigger peak is enhanced above its KS value, and the 
smaller one reduced. This is pole repulsion, and it is felt 
even very far from the strong coupling region. This ef- 
fect is entirely missing from SPA. Next we see that there 
is even a critical value lu c i (d for dark) at which /_ = 
exactly. This means the lower peak disappears entirely, 
and all strength is in the upper peak (Fig. [Sj! From Eqs. 
(dSJ and J2JJ, we find 



An 



SPA 



g(a**)\M 12 \ (A 



0), 



(24) 



where g(a) = 4/tan2a = 16/3 for /f s = 0.1 as is the 
case here. This yields 8.93 eV, whereas the exact result 
is 9.90 eV. 

By increasing u)\ just a little more, we come to the 
position of the avoided crossing uj c (c for crossing) , where 
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Interacting and Kohn-Sham spectra at the critical 
= uJd ~ 9.90 eV. All intensity is in the upper transi- 
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FIG. 6: Interacting and Kohn-Sham spectra for U\ = uj c 
10.61 eV. 
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FIG. 4: Oscillator strengths as function of uj\ . 



9 = tt/2. In fact, Eqs. (JUJ-lHSJ) yield here 

f± = \±(f KS ), (25) 

where (/ KS ) denotes the geometric mean, y/fxf^- ^ n 
our case, this yields /_ = 0.2 and / + = 0.8, respectively, 
giving the lower peak double its KS weight. In FigEl we 
show the spectrum for oj\ — lu c , and observe how much 
it differs from its KS doppelganger. There appears to be 
only one peak, but in fact there are still two, although the 
broadening obscures this. They are very close together. 

The final interesting point is uj e (e for equal), where the 
interacting oscillator strengths equal, i.e., both are 1/2. 
At the equality point, a = 7r/4, and so 8 = tt/2 — 2a KS . 
Again using the high-frequency limit, Eq. 1)2 1[). yields 

An SPA = -4M 12 cot(2a KS ), (26) 



i.e., the same distance above the crossing point, as the 
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FIG. 8: Interacting and Kohn-Sham spectra for u>i = 13 eV. 



amount the point /_ = is below. This yields 12.29 eV, 
whereas the exact number is 11.02 eV (Fig. EJ. 

Finally, in Fig. [SJ we consider uj\ — 13 eV. Now the os- 
cillator strengths have returned (almost) to their KS val- 
ues, but + and — have been reversed. Lastly, we demon- 
strate the dependence of these results on the strength 
of M12 relative to the diagonal elements. We have so 
far presented only the case M\i << Ma. But we have 
argued that it is only the ratio |Mi2|/Af2 SPA that mat- 
ters. Thus increasing Myi does not change the shape of 
the curves (around the turnover point), but only changes 
the scale on which the action takes place. In Fig. ED we 
change Mi 2 to 1 eV and 2.5 eV, and see this occur. Since 
the turnover occurs on a scale of about 4|Mi2|, almost 
the entire region has strong coupling for M\2 = 2.5 eV. 

Lastly, we examine this behavior as a function of 
Mi2- In Figs 1101 and 1111 we repeat the plot of oscilla- 
tor strengths versus u)% for this system, but now with 
Mi2 — 1 eV and M 12 = 2.5 eV, respectively. We see that 
the larger values lead to qualitatively similar behavior, 



but over a broader frequency scale. 

INVERSION 

The above sections present the TDDFT response equa- 
tions in the usual manner. First solve the ground-state 
KS problem, finding occupied and unoccupied levels, 
then calculate matrix elements of the kernel (with some 
functional approximation) , and calculate the true transi- 
tions and oscillator strengths of your system. However, 
we are motivated to gain insight into the excitations, and 
so we ask the reverse question: Given the experimental 
spectrum, what can we learn about the kernel? Inverting 
our equations to solve for 8 yields: 

9 = 2 {a - a KS ) (27) 

Thus, knowledge of the KS oscillator strengths, combined 
only with the experimental oscillator strengths, yields the 
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FIG. 10: Same as Fig. H but for the case M 12 = 1.0 eV. 
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FIG. 11: Same as Fig. but for the case M 12 = 2.5 eV. 



To further extract the matrix elements of the kernel, we 
assume the KS transitions satisfy the same requirement, 
i.e., that the experimental transitions are close to the KS 
ones on the scale of the average transition. This yields: 

M u = (ft- Aft cos 60/4-^/4, 
M 22 = (H + Aftcos0)/4- wi/4, 
M12 = Aft/4sin0. (33) 

These simple expressions give the matrix elements di- 
rectly, once the KS and experimental information is 
known. The mixing angle is determined completely by 
the oscillator strengths, as in Eq. (|27(1 . These expres- 
sions were used to analyze X-ray absorption spectra in 
Ref. [3. 



mixing angle, which measures how strongly the transi- 
tions are mixed! No knowledge of the positions of tran- 
sitions is needed. 

Solving for the diagonal matrix elements we arrive at 



Wu = ^ 2 - (Aft 2 / 2 ) cos ( 

^22 = ft 2 + (Aft 2 /2)cOS( 



(28) 



where ft 2 is the average of ft 2 and Aft 2 is the difference, 
while the off-diagonal matrix element is 



W 12 = (Aft 2 /2) sin( 



(29) 



Again, the experimental positions combined with the 
mixing angle are sufficient to determine the elements of 
the matrix W. The kernel matrix elements themselves 
are then found simply, by using the KS transition fre- 
quencies: 



Mjj = 5k 



and 



Mia = 



Aft 2 sinfl 



(30) 



(31) 



These equations provide an exact way to recover the ma- 
trix elements Wij of the original matrix and therefore 
the matrix-elements My of the kernel K solely from the 
knowledge of the eigenvalues and the angle 6. 

While the above formulas are completely general, in 
practice strong coupling tends to occur between neigh- 
boring transitions. In those cases, the differences between 
the two transition frequencies are often much smaller 
than the transition frequencies themselves. Thus we ex- 
pand in the small parameter Aft/ft, to find 



Wu = ft (ft - Aft cos 6>) , 
W 22 = ft (ft + Aft cos 6>) , 
W 12 = ft Aft sin 8. 



CONCLUSIONS 

To summarize, we have presented the exact formulas 
that arise from a double-pole approximate solution to 
the TDDFT linear response equations. We have shown 
how these reduce to the single-pole approximation when 
the coupling between transitions is weak, and derived the 
leading terms in this expansion, finding results consistent 
with those of Ref. However, with DPA, we can go 

beyond that work, by considering strong coupling. We 
also derive simpler expressions that are valid when the 
transitions are of much higher frequency than the split- 
tings. We illustrated our results, finding (i) that the os- 
cillator strengths can deviate significantly from their KS 
values, even when the coupling is very weak, (ii) that the 
scale to compare the off-diagonal matrix element to is the 
splitting in the single-pole approximation, and (iii) that 
the weaker peak even vanishes at a special value of the 
coupling. 
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